Class 5 Assignment | The Raster Model (v2)

The Raster Model

Concepts & Themes:

This week’s assignment will encompass the following concepts covered in Class 5 lecture & lab:

  • Raster model as primary GIS data model
  • DEM - digital elevation model
  • Continuous vs. discrete rasters
  • Raster surfaces - hillshade, aspect, contours & slope
  • Raster math | reclassification
  • Zonal statistics
  • Vector > Raster and Raster > Vector conversions

Specific Techniques and Tools covered this week will include:

Raster Properties Dialog

Raster Calculator

  • Processing Toolbox Raster Zonal Tools:

Raster Zonal Tools

  • Raster Report for Unique Values:

Raster Report

  • Raster Terrain Analysis Tools:

Raster Terrain Tools

The Raster Model Overview

  • Watch and review this short video overview featuring raster vs vector differences & uses:

Raster Model Overview

Lecture Demo Raster Data Preview:

  • C5 Lecture Demo Overview Sheet

  • C5 Demo Data

  • Tasks:

    1. Previewing Raster Properties
    2. Merging and Clipping Raster Data
    3. Raster Symbolization
    4. Raster Pyramids
    5. Raster Calculator Preview
    6. Discrete Raster Preview (Unique Values Report)

Note: for tasks 2 and 4, the following exercise guide from the QGIS 3.x textbook will be referenced:

Class 5 Readings:

This week’s readings will include 2 sections from the Essentials of Geographic Information Systems textbook; further, 1 supplemental technical readings the raster format within QGIS.

The class 5 quiz will cover only content from the Essentials of Geographic Information Systems textbook as follows:

Class 5 Lab:

  • Class 5 Lab

    • Topics Covered:
      • Raster Calculator
      • Color Ramps
      • Value Tool (QGIS Plugin)
      • Raster Symbolization

Assignment Steps:

  • Assignment Preamble:

In this assignment two (2) map products will be developed - first a thematic map based on gridded (raster) data for global scale population; and second, a topographic mapping based on DEM input data. Both maps incorporate typical raster processing steps integral to early stages of raster analysis.

GPW Overview

HDX Interface - For instance, for Zimbabwe, OCHA produces a dataset containing all 3 admin levels, as it does for many other countries. Choose the zipped shapefile or zipped geodatabase option (both lead to the same page):

HDX Interface

  • On the next page, find the polygon version for all three admin levels. It will likely be towards the top of the list, and the largest file size:

HDX Interface

  • Next, navigate to the zipped file, unzip and connect the .gdb to QGIS:

HDX Interface

  • Note that you now have access to all the admin levels as polygon features necessary of the assignment - admin 0, admin 1, admin 2 and admin 3 geographies. Both admin 0 and admin 2 are necessary for the assignment below:

HDX Interface - Run the raster clipping procedure utilizing admin 0 for your country of choice. This process was demonstrated in the synchronous lecture in a slightly different manner, but same result - vector feature clipping by mask a raster layer:

HDX Interface

Warning: if you receive an invalid polygon error at this juncture, proceed to the section at the bottom of this assignment page titled:

  • Warning - Invalid Geometries in QGIS (workaround Map 1 as needed)

  • Remove the global GPW V4 world pop raster from the layers panel and save the project. Open admin 2 layer from the .gdb if you have not done so yet. You should have 3 files at this juncture - the clipped population raster, the country boundary (admin 0) and the subnational boundaries (admin 2):

QGIS - ADMIN Layers + GPW data

  • Export the admin 2 as a .shp into a project folder you set up for this project. If you do not do this step, you will only be allowed to create virtual fields as the feature will still be inside the .gdb structure.

  • Next, Run Zonal Statistics for the cell value sum per administrative unit. Here the raster cells contain a count of persons per square kilometer. Each of those cells are then summed by Zonal Statistics to the unique geographies of each admin 2 polygon. Access the tool from Processing Toolbox using the search term zonal:

Zonal Stats

  • Populate the tool as follows. Make sure to give a name to the prefix in order to create a new field in the admin 2 feature:

Zonal Parameters

  • The result should appear in the attribute table in similar fashion as follows:

Zonal Results - Sum of Raster Values

  • Next, calculate square kilometers per administrative unit. QGIS uses square meters as default for areal calculations under $area formula in the Field Calculator tool. In order to gain square kilometers, use the following formula in a new field sq.km:

Field Calculator

  • Formula: $area /1000000

  • Next, calculate population density (pop.den) within each administrative unit using the following formula. Make sure to set the output field type as integer (people are represented by whole numbers):

Field Calculator

  • Check your results. The thematic map product will be based on the pop.den field:

Field Calculator - Result

  • Finalize the map product as a thematic map per admin units for population. Use map and legend titling similar to the following:

  • Residents per square kilometer within Zimbabwe Provinces

Note: in the expression above people are named ‘residents’; the unit of measurement is set to per sq. km; and the admin 2 geographies are explicitly named.

  • Utilize skills from class 4 classification methods to choose the best method for your geography. You will want to find a method that best shows the pattern of the data from dense urban administrative units - which tend to be much smaller - to larger rural units.

  • See the map example below to guide your map layout choices. Include a classified choropleth map; corresponding legend, map data source tag and title. Research your country of choice to determine the naming convention for second level administrative units. In the case of Zimbabwe, they are known as districts:

Districts of Zimbabwe

Map Design Example

  • Map II:

  • In this second map, the raster input will be DEM - Digital Elevation Model for wilderness location where a plant species needs both a suitable slope and direction (aspect). When both factors come together, suitability is increased and given a high value. Both slope and aspect can be derived from DEM inputs, which will be done in this second map. The location is the Sandia Mountains on the east side of Albuquerque, New Mexico. To situate this large scale geography, a inset map will be included in the final map layout for the state of New Mexico. To preview, note below the location of the study area within the New Mexico boundary; the attribute name that can be used in the inset map, and a link to a download location for the New Mexico state boundary:

Input DEM situated within New Mexico State

Note: utilize NAME00 for labeling purposes in the final layout map for the map inset. Alternatively, use a text box item from the layout canvas.

Attribute Table View for New Mexico State

  • LINK to New Mexico State Boundary download (use zip file option):

NM State Boundary

  • To start, access/download the assignment PDF and necessary data:

    • Assignment 5, Map II PDF + Data:

  • Load assignment raster data into QGIS:

DEM Loaded to QGIS

  • Navigate to Processing Toolbox, and search Raster terrain analysis:

Processing Toolbox

  • Input for Hillshade:

Hillshade Parameters

  • If you experiment with the Hillshade parameters, your hillshade will be different from the defaults which should look similar to the following:

Hillshade Output

  • Symbolize the DEM, resulting in the following:

Symbolized DEM

  • Switch to the Hillshade, and apply a Blending Mode > Multiply:

Hillshade Properties Dialog

  • The result should be similar to the following. Save the project at this juncture before proceeding:

Result

  • Next, run Slope, and save the result. Make sure to input 35106-B4 DEM, not the hillshade:

Slope Result

  • Next, run Aspect, and save the result. Make sure to input 35106-B4 DEM, not the hillshade or the slope:

Aspect Result

  • Navigate to Processing Toolbox, and begin to develop the reclass table for slope:

Reclassify by Table

Reclassify by Table - Input

  • Populate the table. The book calls for 75.480842590332; however, when run on 3.10 QGIS, a maximum value was found as 76.8519. As a result, the following table has been utilized instead of that suggested by the PDF text. Here the value of 77 is inclusive of the maxiumum value of 76.8519 which will reclass all significant slope values between 55 and the max as 3 :

Fixed Table - Values Populated

  • Output as slopeRC. The result may include a zero 0 value; simply delete that value from unique values as it will not be included in the suitability calculation. Symbolize as follows:

Rendering

  • Check your visual result for slopeRC. Note that you will copy this layer style for the aspect reclass:

Visual Result

  • Save the slope reclass, save the project.

  • Next, reclass by table, but this time with Aspect as input. Populate the table per book example:

Reclass by Table - Aspect

  • Save the aspect reclass result, save the project. Copy/Paste the layer styling from the slope reclass to the aspect reclass per the book example/instructions:

Aspect Layer Symbolized

  • The Final Process Step will be conducted in the Raster Calculator. Here the two suitable surface layers - the suitable slope and the suitable aspect for the plant habitat - will be added together. The result will be values of at least 2 (1 + 1) and no more than 6 (3 + 3). The value 6 will be most suitable, whereas 2 will suitable, but not most suitable:

Raster Calculator

Note: names in the raster calculator above are different than those suggested in the book example. Keep in mind the calculation is as follows:

Slope Reclassification + Aspect Reclassification

  • Check the PlantHabitat.tif results, and symbolize as noted in the PDF for values 2-6. Values increasing from 2 through 6 is the range of suitability where 2 is minimum suitable and 6 is highly suitable.

Plant Habitat Visual Results

Deliverable:

  • Map I and Map II will be delivered to Canvas Assignment 5 upload location as two separate map submissions. Design Map I to include a thematic legend appropriate to the country of choice and its population density pattern (Natural Breaks, Quantile, Equal Interval, ect - your choice). .PDF format submission unless the file size is too large to upload; as alternative, choose .jpg, .png or .TIFF format.

  • Design Map II to include the plant habitat result atop the hillshade raster base. Design a map layout in portrait orientation (right-click inside layout canvas to change paper size and orientation). Utilize the map example to guide your overall design process. Note that this map will feature a custom gradient color ramp for the legend. Utilize the following video guide to develop this additional legend item. Also utilize the New Mexico boundary file to create a titled inset map for the map layout. Export the submission as .TIFF format.

  • Gradient color ramps:

Gradient Color Ramp

Note: for the data source tag, utilize the metadata from the original source for the DEM located HERE

Map II Layout Example

Warning - Invalid Geometries in QGIS (workaround Map 1 as needed)

  • On occassion, QGIS will throw errors due to ‘invalid geometry’. This happens because there is something about the feature - usually a polygon(s) - that QGIS doesn’t see as following all necessary topology rules. Or… and this is the case this week - a country has very small islands that don’t have a valid geometry to conduct the extract by mask process.

QGIS Error Report for Invalid Polygon

Note: Cutline polygon is invalid + point array must contain 0 or >1 elements.

  • One way to repair geometry problems is to simply create a new version of the dataset using a buffering process set to a zero distance. This BLOG POST covers this workaround.

  • However, this process above doesn’t always work, and may not work for this assignment where the offending islands will still be offending in the buffer feature. In this case, a Vector > Raster conversion can be conducted, resulting in a raster layer that can be used as a ‘mask’ when the task is to mask a raster with a vector project area. Here you replace the vector mask with a raster mask utilizing the Raster Calculator.

Step 1:

  • Check your input vector feature. It will likely have multiple polygons, even though the attribute table shows only one feature. This is known as a multi-part geometry and it causes problems when masking rasters if some of the polygons are exceptionally small.

Vector Overlay to Raster Representation

Note: above note the yellow polygons are multiple when overlaid to the population raster. This is a multi-part feature, and some of the features are very small.

Step 2:

  • Convert the country vector to a raster feature. As you do this the input parameters have to be correct as follows:

    • Navigate to the conversion tool:

Rasterize the Vector

  • Input tool parameters as follows:

    • fixed value to burn set to 1
    • input layer = country shp
    • output raster size units = georeferenced units
    • Width/Horizontal resolution = 0.008333333333333329748
    • Height/Vertical resolution = 0.008333333333333329748
    • Output extent = layer extent - choose the country vector
    • Accept all other defaults, including creating a temporary layer.

Rasterize Parameters

  • The result will be a boolean raster with a value of 1 for the country, and no data for everything outside the country.

Boolean Raster Result for Mask

  • Next, open the Raster Calculator. Given that we have a raster of global population, we can overlay that to the new OUTPUTand utilize multiplication. Where 1 exists in the OUTPUT layer, the calculation will result in valid values; everywhere else as no data, outside the frame of analysis:

  • "OUTPUT@1" * "gpw_v4_population_count_adjusted_to_2015_unwpp_country_totals_rev11_2020_30_sec@1"

  • Save the GeoTIFF as 'country_name'_pop

Raster Calculator

Step 3:

  • We can check the results to make sure the new 'country_name'_pop has the same values as the "gpw_v4_population_count_adjusted_to_2015_unwpp_country_totals_rev11_2020_30_sec@1" original layer.

  • First, under symbolization, each layer has to be set to actual values not estimated values.

Set to Actual Values

  • Next, run the Value Tool plugin across the project area; the values between the two layers should be the same, cell to cell. If this condition is met, you can return to the assignment for the next step following the raster extraction step early in the Map I assignment.

Value Tool - Validating Workaround Results